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Abstract 

We carry out a feasibility study for the lattice QCD calculation of the neutron electric dipole 
moment (NEDM) in the presence of the 9 term. We develop the strategy to obtain the nucleon 
EDM from the CP-odd electromagnetic form factor F3 at small 9, in which NEDM is given by 
lim (? 2^Q 0i ? 3(g 2 )/(2m7v) where q is the momentum transfer and is the nucleon mass. We 
first derive a formula which relates F3, a matrix element of the electromagnetic current between 
nucleon states, with vacuum expectation values of nucleons and/or the current. In the expansion 
of 9, the parity-odd part of the nucleon-current-nucleon three-point function contains contributions 
not only from the parity-odd form factors but also from the parity-even form factors multiplied 
by the parity-odd part of the nucleon two-point function, and therefore the latter contribution 
must be subtracted to extract F3. We then perform an explicit lattice calculation employing the 
domain-wall quark action with the RG improved gauge action in quenched QCD at a^ 1 ~ 2 GeV 
on a 16 3 x 32 x 16 lattice. At the quark mass mja = 0.03, corresponding to m n /m p ~ 0.63, 
we accumulate 730 configurations, which allow us to extract the parity-odd part in both two- 
and three-point functions. Employing two different Dirac 7 matrix projections, we show that a 
consistent value for F3 cannot be obtained without the subtraction described above. We obtain 
F 3 (q 2 ~ 0.58GeV 2 )/(2m Ar ) = -0.024(5) e-fm for the neutron and F 3 (q 2 ~ 0.58 GeV 2 )/(2mjv) = 
0.021(6) e-fm for the proton. 

PACS numbers: 11.15.Ha, 11.30.Rd, 12.39.Fc, 12.38.Gc 
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I. INTRODUCTION 



In the strong interaction, one of the most stringent constraints on possible violation of 
parity (P) and time-reversal (T) symmetry comes from the measurement of the electric 
dipole moment (EDM) for neutron (NEDM) and proton (PEDM) d n p . The current upper 
bound is given by 

\d n \ < 6.3 x KT 26 e • cm (90% C.L.) (1) 

for neutronjjj, and 

\d p \ < 5.4 x 10" 24 e • cm (2) 

for protonjij, which are estimated from the results for EDM of mercury atom 199 H g given 
by rf( 199 H g ) < 2.1 x l(T 28 e • cm (95% C.L.) Q 

On the other hand, QCD, which is regarded as the fundamental theory of the strong 
interaction, allows a gauge invariant, renormalizable CP odd operator in the Lagrangian, 
called the 9 term: 

9 f ~ 1 

(3) 

in Euclidean space-time with G^ v the field strength of gluon. Some model estimations 0, 0] 
yield 

|4| ~ 6 x O(10~ 15 ~ 10- 16 )e-cm, (4) 

and this leads to a bound 8<^O(10~ 10 ). Hence 9 must be small or even must vanish in QCD. 

A smallness of 6 in the QCD sector, however, is not protected in the presence of the 
electroweak sector of the standard model, where the quark mass matrix, arising from Yukawa 
couplings to the Higgs field, may be written as 

^(x)M^f(x) + tf{x)M}^f{x), (5) 

where ipi and ipR represent left and right handed quark fields with flavor indices Diag- 
onalizing the mass matrix and making it real, the parameter 9 becomes 

# = #QCD + argdetM, (6) 

where #qcd is the original 9 parameter in QCD. Therefore, either #qcd and argdetM are 
individually small, or the two contributions cancel out to the degree that the experimental 



upper bound on NEDM is satisfied. In either of the two cases, it seems necessary to explain 
why Nature chooses such a small value for 9; this is the "strong CP problem" . One of the 
most attractive explanations proposed so far is the Peccei-Quinn mechanism jf]]. Unfortu- 
nately, the axion, a new particle predicted by this mechanism, has not been experimentally 
observed yet. 

The present model estimations of the hadronic contribution to NEDM, \d n \/9, are maybe 
enough to convince us of the the smallness of #[?], Isj]. However, a theoretically reliable and 
accurate estimation for NEDM will be required to determine the value of 8, if non-zero value 
of NEDM is observed in future experiments. Lattice QCD calculations seem ideal for this 
requirement. Indeed more than 15 years ago the first attempt was made to estimate NEDM 
in a quenched lattice QCD simulation^] by calculating the spin dependent energy difference 
of the nucleon in the presence of the uniform static electric field. In this simulation, the 9 
was converted to the phase of the quark mass term by the chiral rotation, and a non-zero 
value of NEDM was obtained for non-zero 9. Unfortunately this non-zero value of NEDM 
turned outQ] to be a lattice artifact due to the explicit chiral symmetry breaking of the 
Wilson fermion action employed in the simulation; if the 9 term appears in the quark mass, 
an additional disconnected contribution must be included in the quenched calculation of 
Ref.Q. 

After this attempt, there has been no further lattice investigation on this problem for a 
long time. Recently, a new lattice strategy Q| for the extraction of NEDM has been pro- 
posed, and a lattice related investigation on NEDM has been made^|. Also, a preliminary 
lattice QCD result of NEDM estimated from the nucleon electromagnetic three-point func- 
tion has been reported at the lattice 2004 conference jl^]. The value of NEDM divided by 9 
is consistent with zero within statistical errors in this calculation, however. In all attempts, 
an expansion in 9 have been used to avoid the complex action due to the 9 term. 

Since chiral symmetry seems important in the calculation of NEDM^j], it may be prefer- 
able to use lattice fermion formulations having good chiral symmetry, such as the domain- 
wall fermionQ] or the overlap fermionQ]. In these formulations, the bosonic definitions 
of the topological charge agree well with the fermionic definitions. In fact the domain-wall 
fermion was employed in Ref.Q]. 

In this paper we propose a new method to calculate NEDM in lattice QCD. We con- 
sider the strategy of extracting NEDM from the CP-odd part of the electromagnetic form 
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factor of the nucleon, as in the recent lattice calculation [3] as well as in some past model 
calculations [?], Isj]. We first examine the 9 dependence of the three-point (two nucleons 
and the electromagnetic current) function by expanding up to first order in 9 under the 
assumption that 9 has fairly small value. An important finding is that the parity-odd part 
of the three-point function contains a contribution coming from the parity-even form factors 
multiplied by the parity-odd part of the nucleon two-point function. This contribution has 
to be subtracted from the three-point function in order to properly extract the CP-odd form 
factor relevant for NEDM. This point is not considered in Ref. [l3|. We present the formula 
to extract the relevant form factor in terms of the nucleon two- and three-point functions, 
which is valid not only for lattice simulations but also for general cases as long as chiral 
symmetry of quarks is well formulated in the calculations. 

We test our formula in a lattice QCD simulation using the domain- wall fermion in the 
quenched approximation. A renormalization group (RG) improved gauge action is employed 
based on our experience that the residual chiral symmetry breaking due to a finite fifth 
dimensional extent is smaller than in the case for the plaquette gauge action at a similar 
lattice spacing^!]. Since our primary purpose is to check the feasibility of the lattice QCD 
calculation of NEDM based on our formula, we focus on accumulating the statistics to obtain 
a sufficient sampling of topological charge on a 16 3 x 32 x 16 lattice at a~ l ~ 2 GeV with the 
quark mass of rrif = 0.03, corresponding to m^/rrip ~ 0.63. Our numerical results reveal that 
the contribution to the three-point function from the parity-even form factors multiplied by 
the parity-odd part of the nucleon two-point function is really significant. We demonstrate 
the correctness of our formula by showing that two independent ways to extract the relevant 
form factor give consistent results only if the contribution from the parity-even form factors 
is properly subtracted. 

This paper is organized as follows. In Sec. |H] we derive the formula to extract the 
relevant form factor in terms of the two- and three-point functions to first order in the 
expansion in 9. Section ITTT1 contains simulation parameters and technical details. In Sec. II VI 
we present the simulation results for the CP-odd form factors together with the CP-even 
ones. It is demonstrated that we cannot obtain the correct value for the CP-odd form factor 
without subtracting the contribution from the parity-even form factors. We conclude our 
investigation in Sec. |V| 
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II. NUCLEON ELECTRIC DIPOLE MOMENT FROM CORRELATION FUNC- 
TIONS 



A. Form factors and EDM 



We consider the electromagnetic form factors of the nucleon, defined by 



(N(p,s)\JnN(p',s')) = u(p,s) 



2m 



N 



u(p', s'^ 



(7) 



where q = p — p' is the momentum transfer, \N(p, s)) is the on-shell nucleon state with 
momentum p, energy p = m 2 N + p 2 and helicity s. The electromagnetic current for 



quarks, J™, is given by 



(8) 



/ 

where ipf is the quark field with flavor / and the electric charge e/. 

In the small momentum transfer limit, q — > 0, the above form factor is described by the 
following effective interaction: 

"*s(0) 



N(x) 



2m 



N 



N(x), 



(9) 



where A^(x) is the U(l) electromagnetic field. We can define the electric dipole moment 

as 



\d 



A' 



lim 



F 3 (q 2 ) F 3 (0) 



(10) 



q 2 ~*o 2m^ 2m jy 

So far we do not specify CP properties of the system and therefore dependences are omitted 
in the above expressions. 



B. Extraction of form factors from correlation functions 

In the lattice calculation we must extract the form factors from the three-point correlation 
functions such as 

G e NJ ^ N (q,t,r) = (e\N(p,t)J™(q,T)N$,0)\e), (11) 

where N(p, t) and N(p, t) are the interpolating fields of the nucleon at time t, which contain 
annihilation and creation operators with momentum p. Here \9) is the 9 vacuum, which may 
be defined as 

\9) = e l9Q5 |0), Q 5 = J d 3 xK (x,t), (12) 
6 



using the conserved (but gauge non-invariant) chiral current K^. Throughout this paper, 
we do not use the explicit representation oi\9). Inserting the complete set of states between 
the interpolating fields and the current, we obtain 

G e NJ ^ N (Q,t,r) = e-s^e-*** 

x J>|iV|iVV, s))(N e (p, s)\J™\N e (p', s'))(N e (p>, s')\N\9) + ■ • • , (13) 

s,s' 

where E N e = ^Jp 2 + m 2 N0 , E' N6 = \jp' + m 2 N0 , and the dots represent exponentially sup- 
pressed contributions, which are assumed to be negligible. Hereafter we represent explicitly 
superscript of 9 on 9 vacuum effects. Since the interpolating fields N and N can be expanded 

as 

N(p, 0) = Z% a N° (p, s)u%(p, s) H , 

s 

N(p, 0) = (Z e N y <J N M s)u 6 N fa *) + ••-, (14) 

s 

we have 

(6\N\N e (p,s)} = Z e N u e N (p,s), 
(N e (p',s')\N\9) = (Z e N yu%(p>,s>), (15) 

since 

\N e (p,s)) = a^ Ne (p,s)\6) (16) 
(N e (p,s)\ = (9\a N e(p,s). (17) 

In eq. (fTlj) the dots represent contributions from all other possible states including the 
negative energy state of the nucleon, the negative parity states or other excited states. 
Except a special case considered in the appendix, we will not consider these states in this 
paper. The spinors u e N (p,s) and u e N (p,s) are on-shell nucleon wave functions which satisfy 
the Dirac equation: 

(27 ■ p + m N ee- lfN{eh5 ) u%(p, s) = u%(p, s) (ij-p + m N ee- ifN{0h5 ) = 0. (18) 

Since CP is broken in the 9 vacuum, a CP non-invariant phase factor e l ^ N ^ 15 can appear in 
the mass term. This means that /tv(6') must be odd in 9, while m N and Z N are even in 9. 
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Therefore we have 

m N e = m N + 0(8 2 ), Z e N = Z N + 0{8 2 ) 

fs(0) = fho + oie 3 ) (19) 

for small 9. The bi-spinor projection satisfies 

x-e/- x -i 1 -p + m Ne e i M e ^ 
yu N {p,s)u N {p,s) = — (20) 

s 

with the normalization that u e N (p, s)u e N (p, s') = — — — _ ^ g ,_ 

The form factor, which we wish to extract from the three-point function, is written as 

(N°{p, s)\J™\N e (p>, s')} = u%(p, s)W^q)u%(p', s>) (21) 

It satisfies 

4r(P, s)q»Wl(q)u e N {p\ s') = (22) 

from the current conservation d^J^ M = 0. Taking account of the parity properties, the form 
factors in general should have the form that 

= 9{0 2 ) W;™(q) + id h{8 2 ) W° dd (q), (23) 
where we take g(0) = h(0) = 1 for their normalizations and 

W^(q) = ^(q 2 ) + (24) 

Wf d {q) = ^\<V75 + F A (q 2 ) (q^ ■ q - ^q 2 ) 7s (25) 

with a^ v = iljn, 7^/2. Here Fi^{q 2 ) are the electromagnetic form factors and FA(q 2 ) is called 
"anapole form factor". As mentioned in Sec. Ill Al F 3 (0)/(2m7v) gives the nucleon electric 
dipole moment. For small 8, we can expand g(6 2 ) = 1 + 0{8 2 ) and h(8 2 ) = 1 + 0(8 2 ). 
Using all the necessary expressions and expanding in terms of 8, we obtain 

-27 -p + m N {l +ifjf8j 5 ) 



G%jJq,t,T) = \Z N \ 2 e- E ^e- E '^- 



2E 



N 



x [w r(q) + Uiwrm dULt^Jim + 0(S »). (26) 

On the other hand, G e NJ ^ N (q, t, r) can be evaluated by the path integral as 

G% JtlN (Q,t,r) = (N(p,t)J^(q,r)N(p',0)e^) 

= G NJ ^ N (q, t, r) + i9G Q NJuN {q, t, r) + 0{8 2 ) (27) 



with 

G NJliN (q,t,r) = (N(p,t)J™(q,r)N(p>,0)), (28) 
G% JfiN (q,t,r) = (N(p,t)J™(q,T)N$,0)Q), (29) 

where N, N, and J™ are c-number operators in the path-integral, albeit using the same 
notation as the quantum operators, and Q is the topological charge, defined by 

Q = ^S f d 4 xG, u (x)G, u (x) (30) 



32tt 2 „ 

in the continuum theory. The symbol (O) denotes the path-integral average of an operator 
O in QCD with 9 = 0, which is given by 

(O) = J Vi)Vi)VA^ e S G (A)+i>D{A)^ 0) ( 31 ) 

where Sc(A) and 2jjD(A)ip are the QCD gauge and quark actions, respectively. 
By comparing eqs.(J2EJ) and ()27j1 we obtain 

t,r) = \Z N ?e-^>-%-°»< ■' + m " w r{q) =*? /', + m » , (32) 



G^iqXT) = \Z N \ 2 e- E ^e-^ 0J 7 W° dd (q) 



-ry • p + m N dd -27 • p + 



2-^v M 2E'n 

fk^N TT/even/ x ~»7 ' + m N 

2E N 75 /1 W 2£^ 

-7-p + -,^ rn(g) ^ 7 l (33) 



2.EW m 

These formulae are one of the main results of this paper. They show that the parity-even 
form factor W° vcn multiplied with the parity-odd term /jy of the two-point function has to 
be removed from the three-point function G^j^j at order 9, in order to obtain W° dd (q). We 
will see in Sec. II VI how significant this subtraction is. 



C. Extraction of fjj from the nucleon two-point function 

In order to determine /jy we consider the two-point function of the nucleon given by 
G% N {p,t) = (S\N(p,t)N(p,0)\9) 

2E N e 

= (N(p,t)N(p,0)e id ^). (35) 
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By comparing the operator expression ()34|) with the path integral expression (J35|) at each 
order of 6, we obtain 

G NN (p,t) = (N(p,t)N(p,0)) = \Z N \ 2 e- E ^ ~ %1 ' P _ +mN (36) 

at 8° and 

G Q NN (p,t) = (N(p,t)N(p,0)Q) = |^| 2 e-^4^75 (37) 

at 9 l . These formulae tell us that can be determined numerically from the nucleon 
propagators through the appropriate insertion of Q. 



III. SIMULATION DETAILS 
A. Fermion and gauge actions 

For the quark field we employ the domain- wall fermion actional : 

S F = - 2^ ^( x ' s)D DW (x, s; x', s')ip(x ', s') + m f q(x)q(x), (38) 

x,s,x',s' X 

D DW (x, s; x's') = D\x, x')5 s , s , + 6 x ^D 5 (s, s') + (M — 5)<5 X|X /<J a>8 /, (39) 
D 4 (x, x') = £ ~ [(1 - 7^0*) <W + (! + 7^(z')<W+/J , (40) 

A* 

i(l-75)V (3 = 1) 

D 5 (s, a') = < |(1 - k)5 s+ i, s , + |(1 + 7s)<W (K s < 7V S ) , (41) 

|(l + 75)^-i )S ' = 

v 

where x, x' are four dimensional space-time coordinates, s, s' are coordinates in the fifth 
dimension of length N s , and M is the domain- wall height. Here q, q are four dimensional 
physical quark fields given by 

q(x) = ^(l- l5 )iP(x,l) + ^(l + l5 )^(x,N s ), (42) 

q(x) = ^( z ,JV.)i(l-75)+^(x,l)i(l + 75) (43) 

and m/ is the bare quark mass. This action has exact chiral symmetry at N s — > oo, and 
therefore it is expected that the fermion fields respond properly to the topological charge of 
the gauge fields. 
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For the gauge field we employ the following improved action: 



L plaq rect ) 



where the first term represents the plaquette and the second the 1x2 rectangle loop. In this 
paper we take en = 3.648 and c\ = —0.331, proposed by Iwasaki using the renormalization 
group analysis (l8j. The action with this parameter choice is called the RG improved gauge 
action. 

B. Topological charge 

In continuum QCD, the topological charge density is defined by 



<l( x ) = TTTr^^Tr G^ v (x)G a p(x) (45) 

OZ7T z 



and the topological charge is given by the integration of q(x) over space-time: 

Q = (46) 

In our lattice calculation we define the topological charge by the improved bosonic 
formal, which is given by 

Qimp = Y [^QplaqO) + CiQ Iect (x)] ( 47 ) 
x 

I p p 2 R R 

Qplaq = 22^2 e livap r R(F IJiV F a { 3 ), Q rcct = -^j^E ^a^Tl (F^F a/3 ) (48) 

where F^ v and Fj^ u are field strengths constructed from the plaquette clover and the rectan- 
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gular clovers as depicted in Fig. ^ Their explicit expressions are given by 



Im 



— ^ tr \Ufa)U v {x + fl)Ufa + v)Ufa) 



X,IM,U 



+ U v {x)Ufa - pL + v)Ufa - flU^x - pi) 

+ Ufa - p)Ufa — ft - v)Ufa - pi - 0)U u (x - u) 

+ Ufa - v)Ufa - i>)U u {x + pi- v)Ufa) 
Im 



(49) 



x — v) 



F Hi*) = y ^ tr [[/.(^^(x + + /i + ,))(/J(x + 2/))f/](x + i))f/f (x) 

+ U v {x)U v {x + u)Ufa - fi + 2v)Ufa - fi + i>)Ufa - p,)Ufa - ft) 
+ Ufa - p)Ufa -pi- i>)Ufa -pi- 20)U ll {x -pi- 2v)U v {x - 2v)U v { 
+ Ufa- u)Ufa - 2u)U„{x - 2u)U v {x + pi- 2u)U v {x + pi - u)Ufa) 
+ tr Ufa)Ufa + pi)U v [x + 2pL)Ufa + v + pL)Ufa + v)Ufa) 
+ C^(z)yJ(z - pi + v)Ufa -2pi + v)Ufa - 2pL)Ufa - 2p)Ufa - pi) 
+ Ufa - pi)Ufa - 2pt)Ufa -2pi- u)U^{x -2pi- u)U^x -pi- v)U v {x - u) 
+ Ufa - 0)U^(x - v)Ufa -u + pi)U u (x + 2pi- v)Ufa + pi)Ufa) j. (50) 

In order to remove 0{a 2 ) discretization errors in Qim P at the tree level, we take c\ = —1/12 
and c = 1 - 8ci = 5/3I20J. 

C. Nucleon three-point functions and form factors 



To compute nucleon three-point functions, we construct the amplitudes for the six di- 
agrams in FigEl each of which consists of two quark propagators and one double quark 
propagator with an insertion of the electromagnetic current J™: 



J™{x) = e u V;(x) + e% d (x) + e%*(*) 
with Vfa) (i = it, d, s) the conserved vector current of flavor i defined by 



Vfa) 



E 



JuK x > s ) 



1 r T 



ip l (x, s)(l - 7 M )£^(x)^*(x + pt, s) 



+ft(x + pi, s)(l + ^)Ufa)^(x, s) 



(51) 



(52) 



(53) 
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Since this current is conserved, the renormalization factor satisfies Zy = 1. To insert V?, we 
use the ordinary source method. 

We do not consider the disconnected diagrams depicted in Fig. El since this contribution 
vanishes in the SU(3) flavor symmetric limit due to the fact that e u + e d + e s = 0. In future 
investigations, however, it is necessary to check if this contribution is indeed small or not 
for m U)d ^ m s . 



D. Lattice parameters 

We generate quenched gauge configurations on a 16 3 x 32 lattice with the RG improved 
gauge action at (3 = 2.6. Gauge configurations separated by 200 sweeps are used for mea- 
surements after 2000 sweeps for thermalization, where one sweep consists of one heat-bath 
update followed by four over-relaxation steps. The quark propagator is calculated using the 
domain-wall fermion with N s = 16 and M = 1.8, where N s is a size of the fifth dimension 
and M is the domain- wall height. These parameters are the same as those employed in 



Ref. 2jJ, in which the lattice spacing was determined to be or 1 = 1.875(56)GeV from the p 
meson mass. The residual quark mass m ros , which represents the explicit chiral symmetry 
breaking effect at finite N s , is found to be m res ~ 4MeV. 

We take the bare quark mass rrifa = 0.03 corresponding to the ratio of pseudoscalar 
to vector meson masses mps/my = 0.629(8) [2^. To reduce the effects of unwanted ex- 
cited hadron states, we calculate the quark propagator employing an exponentially smeared 
source, given by Ae~ Bra with A = 1.28 and B = 0.40, where ra is the distance from the 
origin at (7, 7, 7) in the lattice unit. We take the periodic boundary condition for both 
gauge and fermion fields in order to minimize possible violation of the equivalence between 
the bosonic and the fermionic definitions of topological charge. 

In the calculation of three-point functions, the electromagnetic current is inserted at the 
fixed time slice r = 6. The nucleon sink is varied as a function of t with the nucleon source 
fixed at t = 1. We move the final state nucleon with one of the three non-zero spatial 
momenta p = (7r/8,0,0), (0,7r/8,0), (0,0,7r/8), while the initial state nucleon is placed at 
rest p' = (0,0,0). 

We accumulate 730 configurations for measurements. Errors are estimated by the single 
elimination jackknife procedure for all measured quantities. 
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IV. NUMERICAL RESULTS 



A. Topological charge 

In order to reduce ultra-violet fluctuations of gauge configurations in the calculation 
of the topological charge, we employ the cooling method using the action in eq.([44[) with 
C\ = 3.648 and c = —0.331. In Fig. HJwe show both Qpiaq and Q irap as a function of the 
cooling steps, where Q p i aq is the naive plaquette definition of the topological charge and 
Q imp is the improved one in eq.(|4T|). It is clearly observed that Qi mp approach integer values 
much faster than Q p i aq as the cooling step increases. The values after 20 cooling steps are 
also closer to integer values. The deviation from integer is within a few percent for Qi mp 
after 20 cooling steps. As the topological charge, we employ Q imv measured after 20 cooling 
steps on each configuration. 

In Fig. |S] we show the time history of the topological charge as a function of the gauge 
configuration separated by 200 sweeps, where we do not recognize any long autocorrelation. 
In addition, the histogram of the topological charge in Fig. |B] exhibits nearly a Gaussian 
distribution. These observations suggest that the sampling of the topological charge is 
sufficiently good in our calculation. 

B. Nucleon propagator 

1. 0(9°) contribution 

From the analysis in Sec. Ill CI the 0(8°) contribution of the nucleon propagator should 
behave as 



where T4 = (1 +74)/2 and the dots represent the contribution from the excited states. Here 
\p\ = 1 is a shorthand notation of p = (ir/8)n with \n\ = 1. 

The nucleon propagators from the smeared source and the point sink with \p\ = and 
\p\ = 1 are plotted in Fig. [7J and the corresponding effective masses are given in FigJHJ 




(bl = o) 
(bl = i) 



(54) 
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Choosing 9 < t < 13 for the fitting range, we obtain 



E^a 




0.7114(43) {\p\ 
0.8068(77) 



2162(87) (\p\ 
1453(104) (\p\ 



0) 
1) 



0) 
1) 



(55) 



Note that E^a = 0.8068(77) at \p\ — 1 agrees with the relativistic dispersion relation 
y/(m N a) 2 + (vr/8) 2 = 0.813(4) within the error. 

2. 0{9 l ) contribution and the determination of f N 

In Sec. Ill CI we have shown that the parity-odd part of the nucleon propagator should 
have the following form at the order of 9 1 : 



In Fig El we plot the parity-odd part with the insertion of Q, while the parity-even part 
without Q, tr [Gnn(p, £)y], is given in Fig [IT]] for comparison. The latter vanishes at all 
time slices as expected. On the other hand, the former is non-zero, showing an exponential 
fall-off at large t. After averaging over the forward and backward propagators in time, we fit 
the party-odd part at 9 < t < 14, by the form £e~ Erft with En fixed to the value obtained 
from the previous fit of the 0(6°) contribution. The fit, represented by the solid line in 
Fig. El is reasonably good, showing that the decay rate of the parity-odd part is consistent 
with that of the parity-even part for both \p\ = and 1. Note however that the small 
deviation indicates that the present statistics are still insufficient at large t for \p\ — 1. 

Another confirmation of consistency is found in the effective mass plot of Fig. ^2 where 
the parity-odd part of the nucleon propagator with Q and the parity-even part without Q 
are compared. 

With the aid of En and Zn obtained previously, £ is converted to f N as 



tr [G NN (0,t)f] 
tr [G Q NN (p,t)^] 



7 \ 2 f X P~ 



(56) 



Zn(p)\ 




(57) 




(58) 
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Despite the deviation at large t for \p\ = 1, the two estimates of f N agree within the errors. 
In our analysis we will use f N obtained from \p\ — 0. 

We note that the exponential fall-off of G NN is very sensitive to the distribution of 
the topological charge. One could use this quantity to judge whether the sampling of the 
topological charge is sufficiently good: if the exponential fall-off of G^ N does not agree with 
that of Gnn, the sampling is not satisfactory. 



C. Form factors 



1 . Extraction from three-point functions 



In order to remove the exponential fall-off due to the nucleon propagation in the three- 
point functions, we divide them by the two-point function Gnn(p, t) with the smeared source 
at t — 1 and the point sink at t: 

with 

lg ' T ' } [tr[r 4 G^(0,t)]tr[r 4 G^(0,r)]tr[r 4 G^(g,t-r)]J ' 1 ' 

where the factor R is to remove the interpolating field dependent normalization factors. 
Here we introduce a new two-point function G N p N (q,t) that has the point source at t = 1 
and the point sink at t, and G N ^(q,t — r) represents that the source point is started at 
t = t. In Fig. [T21 we plot the t dependence of R(q;r,t) with r = 6 and \q\ = 1. As we 
expect, the correction factor R is independent of t in the large t region. 

The form factors are extracted from the ratios defined in eqs. ()59|) and (JoDj) by apply- 
ing appropriate projections of gamma matrices. From the 0(6°) contribution H\(q]t, r) of 
eq. (|59|) . F\ and F 2 can be extracted as 



E N \ Am~ N 



tr[n^;t,r)r 4 ] = ~" ' I Ftf) + F 2 (q^-^ , (62) 



tr [n*(q; t, r)zT 47 57i] = ^ J 1 {HQ 2 ) + F 2 (q 2 )) , (63) 
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where we use eq. 
form factors as 



The F\ and F 2 are related to the electric (G e ) and the magnetic (G m ) 



G e (q 2 



.2} Q 



F 1 (q 2 ) + F 2 (q 2 



G m {q 2 ) = F 1 (q 2 ) + F 2 (q 2 ). 



(64) 
(65) 



These quantities become the electric charge and the magnetic moment at q 2 — > 0. Therefore 
G c (0) = 0(1) and G m (0) = -1.91 (2.79) [22J for the neutron(proton). 

We can extract the form factor F^(q 2 ), which is relevant to the electric dipole moment, 
in two different ways from H1(q; t, r) of eq. (pH|) : 



tr [n°(g,t,r)r 475 ] 



tr [n°(g,t,r)zr4757i] 



Fz(q 2 



+ 



2E N m N 
E N + m N 



2E N 
En + m N 
2E N m N 



FM 2 ) + 



AmiyEN 



-F,{q 2 



q t F 3 



q, 



-Friq 2 



qi{E N + 3m N ) 



F 2 {q 2 



N- 



(66) 



(67) 



2E N Am N E N 
where we use eq. fl33|) . It will be clear that the consistency between the two different extrac- 
tions of Fs(q 2 ) provides a crucial test for the validity of actual calculations. Note that the 
anapole form factors which appears in W° (q) is canceled out after applying the trace of 
spinor with these projections. 



2. Parity- even part 

We first consider the parity-even part of the form factor. In Figs. EI] and El we show 
the magnetic and electric form factors obtained from eqs. ()62j) and (|63jl . We observe good 
plateaux except for the electric form factor of the neutron. Employing a constant fit with 
10 < t < 14 for G m and 12 < t < 15 for G e we obtain 

f -0.591(37) (neutron) 

G m (q 2 ) = { (68) 

[ 0.952(60) (proton) 
G e (q 2 ) = 0.502(33) (proton) (69) 

at q 2 = 0.58 GeV 2 . Here the value for G™ is not quoted since it is difficult to find a plateau 
for the fit. Note also that the value of G e slightly depends on the fitting range, as could be 
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seen from the figure. These values are compared with the previous lattice result |23| obtained 
at a similar nucleon mass and lattice spacing using the nonperturbatively 0(a) improved 
Wilson quark and the plaquette gauge actions: G m (q 2 ) = —0.73(4) for the neutron, and 
G m (q 2 ) = 1.17(6) and G e (q 2 ) = 0.451(13) for the proton. 

Experimentally the proton electromagnetic form factors at small q 2 have been measured 

as | 

G p Jf* p \q 2 = 0.58GeV 2 ) = 0.848 ± 0.011, G p e {cxp \q 2 = 0.58GeV 2 ) = 0.294 ± 0.008. (70) 
Recent experimental measurements for neutron form factors cover a wide range of q 2 from 



251 ]. For example, the 



very small q 2 ~ 0.3 GeV 2 to more than 1 GeV 2 with high precision 
neutron magnetic form factor at q 2 = 0.6 GeV 2 is measured as 2(| 

G" (cxp) (<? 2 = 0.6GeV 2 ) = -0.568 ± 0.007(stat.) ± 0.015(syst.), (71) 

and the electric form factor at q 2 = 0.5 GeV 2 isj^ 

G n(cx P )^2 = .5GeV 2 ) = 0.0463 ± 0.0062(stat.) ± 0.0034(syst.). (72) 

The fact that our results are not so much different from these experimental values may 
suggest that systematic errors involved in our calculations, such as the large quark mass, 
the small lattice size, the non-zero lattice spacing and the quenching, are not so large on 
these quantities. 

The form factors Fx and F 2 are obtained from G m and G e using eqs. (|64j) and (J65|) . We 
show the time dependence of Fx and F 2 in Figs.E3andE3 We observe plateaus at t > 13 for 
F 2 n ' p and Ff but not for F™, which behaves badly as does. We think that this behavior is 
caused by the contribution from other excited states, as pointed out in Ref.j3]. Therefore 
we have tried to include the contributions of the parity-odd nucleon state AL in our analysis, 
using the formulae given in Appendix [XJ Although larger errors prevent us from extracting 
the reliable values for G™ and F™, we find a much better plateau for G" as shown in Fig. [T7| 
and for F" in Fig. [THJ 

Applying a constant fit to Ff and F% ,p with 13 < t < 16, we obtain 

F p (q 2 ) = 0.515(37), F^q 2 ) = -0.560(40), F p (q 2 ) = 0.399(37), (73) 
which are compared with Ff{q 2 ) = 0.499(13), F p (q 2 ) = 0.68(6) fron 
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3. Parity- odd part 



We now consider the main target of our calculation: the parity-odd part of the form 
factors. Substituting G®j N (q; t, r), Fi(q 2 ), F 2 (q 2 ) and /jy into eqs.([6bj) and (J57|) . we ex- 
tract the form factor F 3 (q 2 ) in two different ways. We first plot tr [Il^g, t, r)r47 5 ] and 
tr [Il^g, t, r)T4757j] in Figs. for the nucleon and in^Ulfor the proton, together with the 
second terms in the right hand side of eqs. (J5fiJl and (J57|) . composed of Fi, F 2 and f^. Com- 
paring the two extractions, we find that tr [n^g, t, r)r 4 7 5 ] and tr [IT^g, t, r)r 4 7 5 7j] shown 
by open symbols do not agree with each other. In particular, for the case of the proton, 
tr [n^g, t, r)r 4 7 5 ] is very large, while tr [Il^g, t, r)r 4 7 5 7j] is almost zero. 

We note that tr [Il^g, t, r)r 4 7 5 ] and tr [n^g, t, r)r 4 7 5 7j] are almost zero for the neutron 
is consistent with the previous result ^J, in which the contribution of the mixing term was 
not considered. 

Figure l?T1 shows the result after subtracting the second terms from tr [Il^g, t, 1^475] 
(circles) and tr [Il^g, t, r)r 4 7 5 7j] (triangles) for the neutron (top panel) and the proton 
(bottom panel). It is gratifying to observe that the two ways of projections give values 
of F 3 consistent with each other. In particular, the significant difference between the two 
projections for the proton almost vanishes after removing the second terms. 

Making a constant fit to the result of eq. ()67|) over 12 < t < 15, we obtain 

1 [ -0.024(5) e ■ fm (neutron) 

F 3 (q 2 ~0.58GeV 2 ) = <^ K J 1 ; (74) 



2m N 0.021(6) e-fm (proton) 



with a^ 1 = 1.875 (56) GeV(2l|. Note that the opposite sign between neutron and proton 
is consistent with the estimate of Ref.jsJ, where the baryon EDM is proportional to its 
magnetic moment. It should be remembered that our values are obtained at a non-zero 
q 2 . With the same caution in mind, our result may be compared with other estimates of 
F 3 (0)/2m7v: 3.6 x 10" 3 e-fmQ, 7.5 ± 3.2 x 10- 3 e-fm|s| and 3.9 x lO^e-fmJ 



V. CONCLUSION AND DISCUSSION 



In this paper we have carried out a feasibility study for the lattice QCD calculation of 
the neutron electric dipole moment in the presence of the 8 term. We took the strategy 
to extract the nucleon EDM from the nucleon electromagnetic form factor at small 9. We 
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explained how one could extract the parity-odd part of the form factors, which becomes 
EDM at zero momentum transfer, from the relevant correlation functions for small 9. We 
have pointed out that the contribution of the parity-even form factors multiplied by the 
parity-odd part of the nucleon two-point function has to be subtracted from the parity odd 
three-point function. We have derived the formula to carry out this subtraction, which is 
one of the main results. 

In the second half, we applied our formula to an actual lattice QCD calculation, employing 
the domain-wall quark action with the RG improved gauge action in quenched QCD. We 
used a 16 3 x 32 x 16 lattice at a -1 = 1.875(56)GeV with the quark mass of rrif = 0.03, which 
corresponds to m v /m p = 0.629(8). In order to obtain a sufficient sampling of topological 
charges, we accumulated 730 configurations, which enabled us to extract the parity-odd 
part both in the nucleon two- and three-point functions. We have shown that two different 
7 matrix projections yield values for F 3 consistent with each other, if the subtractions 
mentioned above are made. 

Our lattice calculation gives F 3 (q 2 ~ 0.58GeV 2 )/(2m7v) = -0.024(5) e ■ fm for neutron, 
Fs(q 2 ~ 0.58 GeV 2 )/(2m7v) =0.021(6) e-fm for proton at m^a = 0.03. In order to obtain the 
physical value of NEDM, we have to perform various extrapolations: q 2 — > extrapolation, 
the chiral extrapolation (rrif — > m p h ys ) and the continuum extrapolation (a — > 0). In 
addition, we have to remove the quenching error in NEDM, since NEDM vanishes at zero 
quark mass in QCD but it does not in quenched QCD. While it would be possible to remove 
each of these systematic errors step by step using the method proposed in this paper, we 
are now looking for better alternatives with which we can remove them more easily. 
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APPENDIX A: JV. TO iV FORM FACTOR 



If we include the electromagnetic transistion from N + to AL, eqs. f!62|) and ()63|) are mod- 
ified as 

,2 



where 



tr[n!(^,r)r 4 ] = EN + mN [ Ftf) + F 2 (q 



E 



N 



„2\ 1 

4m 2 



N 



+ 



+ 



E N _(m N + m N _y 3 
q 2 (E N _ - m N _)(l - A) ^jv+^at 



E 



N- 



A+_(t,r) 



tr[n j e (g;t,r)zr 4 7 5 7 i ] = e ijk -^- 



7$> 



QkQ T-iN + ->N-r~2 
+ — — r 



E 



(Fife 2 ) + F 2 {q 2 )) 
qk{m N -m N _) N 

7 \~ 3 

E N _(m N + m N _) 

A+_(*,r), 



(m N + m N _)(E N _ -m N ) 
A = — - — , A+_(t, T) 



Z, 



± e -(^iv_-^iv)(*-^) 



and m^v_ is the parity-odd nucleon mass and F^ + ~* N ~ (q 2 ) and F^ + ^"~ (q 2 ) are the N + to 
AL transition form factors with q 2 = —{E^_ — m^v) 2 + (f ■ (p) is the parity-odd nucleon 
amplitude which can be determined from the AL propagator. In order to obtain these form 



(Al) 



(A2) 



(A3) 
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factors we have to consider additional projection formula as follows: 
tr[n|(g;t,r)zr 4Tj ] = ^ Utf) + ^-F 2 (q 2 ] 



+ 



+ 



qj(E N _ + m N _] 
En Amu + m N _ 



iN+—>N- / ~2\ 



m 2 



E 

Qi 



A + _(*,r) 



(A4) 



tr[n j e (g;t,r)zr 4 ] = ^-Ftf 

&N 



qi{E N -m N ) 

2i71nEn 



F 2 (q 2 ) 



+ 



qi(EN^-m N ) z - 1 jv+^jv_/~2n 
^3 [Q ) 



(m N + m N _) E N. 

?N 



qiq 2 (E N _ -m N + {-E N _ + m N _)A) z? jv + - + jv_/~ 2 ^ 



tr[n^(g;t,r)r 47j ] 



En— (E n _ - m N ) 
-E N + m^v 



A+_(*,r) 



(A5) 



E 



N 



(E N - m N ) 2 Sij - (fSjj + g,gj c 
2E N m N 

{E N _ + m N _)(E N _ - m N )8ij - + N+ ^N- r , 



+ 



E 



(m N + m N _)E N _ 
-{E N _ +m N _)5ij + 



Aq^j 



E N _ - m N 



,7V+->7V_,~2^ 



A + _(t,r). 

(A6) 



Employing four equations out of the above five formula we can determine Fx, F 2 , p^+^ N 
The form factors F\ and F 2 plotted in Figs. El and are obtained with the 



and F 



N + ^N- 



use of eqs.tjXU), (jMJ) and $Kty(i = j and i ^ j). 

In Fig. 123 we show the effective mass plot for the parity-odd nucleon. As for |p| = 1 
we have subtracted the N + contribution from tr [(jjvjv^jP)^] with F 4 = i^r 1 , in order to 



2 ' 



obtain the energy En_ for the parity-odd state. The global fit with 4 < t < 9 for rriN_ and 
6 < t < 11 for En„ gives 



m N _a = 1.017(8), E N _ a = 1.056(18) 



(A7) 



where En_cl is close to the value expected from the relativistic dispersion relation 



^{m N _a) 2 + (vr/8) 2 = 1.090(8). 
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FIG. 1: Plaquette clover (left) and rectangular clovers (middle and right). 
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FIG. 2: Six diagrams contributed to the neutron three-point function. Cross symbol represents 
the current insertion. In (a), (6) EM current is connected with u quark line, and in (c) ~ (/) with 
d quark line. 



26 



11 



u, d, s 




11 



FIG. 3: An example for the disconnected contribution to three-point function. The loop contains 
u, d, s quarks 
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FIG. 4: Topological charge as a function of the number of cooling steps. Open symbols represent 
the naive plaquette definition and solid symbols the improved definition. 
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FIG. 5: Time history of the topological charge. Gauge configurations are separated by 200 sweeps. 
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FIG. 6: Histogram of the topological charge. 
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FIG. 7: tr [Gnn(p, tyT^], 0(0°) contribution of the nucleon propagator with T4 projection for 
\p\ = (top) and \p\ = 1 (bottom). Solid lines denote the fitting results. The parity-odd state 
dominates the backward propagation in t denoted by triangles. 
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FIG. 8: Effective mass plot for the nucleon with \p\ = (top) and \p\ = 1 (bottom). Solid(dotted) 
lines denote the central values(errors) of the global fit. 
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FIG. 9: tr Gj rN (p,t)-^- , parity odd part of the nucleon propagator with Q for eq. (|56j) |p| = 
(top) and eq.(|57[) |p| = 1 (bottom). The sign of propagator is negative. Solid lines denote the 
fitting results. 
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FIG. 10: tr [Gnn{p, > parity odd part of the nucleon propagator without Q for eq. Qftfijl \p\ = 
(top) and eq.(JS7J) \p\ = 1 (bottom). 
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FIG. 11: Effective mass plot for the parity odd part of the nucleon propagator with Q (filled) for 
ea.ea. ()56|) \p\ = (top) and ea.(|57|l \p\ = 1 (bottom) after averaged with backward propagator in 
time direction. The parity even part without Q (open), which are already shown in Fig. IHl are also 
plotted for comparison. 
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FIG. 12: Correction factor R(q; r, t) as a function of t with r = 6 and \q\ = 1. 
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FIG. 13: Magnetic form factor for neutron (top) and proton (bottom). 
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FIG. 14: Electric form factor for neutron (top) and proton (bottom). 
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FIG. 15: Form factor Fi(q 2 ) for neutron (top) and proton (bottom). Horizontal lines represent the 
constant fit. 
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FIG. 16: Form factor F2(q 2 ) for neutron (top) and proton (bottom). Horizontal lines represent the 
constant fit. 
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FIG. 17: Solid circles represent the electric form factor G™(q 2 ) for neutron in different analysis 
where the effect of iV_ is taken into account, together with the previous result(open) in Fig|14l 
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FIG. 18: Solid circles represent the form factor F™(q 2 ) for neutron in different analysis where the 
effect of N_ is taken into account, together with the previous result(open) in Fig |15l 
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FIG. 19: Parity-odd part of the form factors for the neutron. Top: tr [n^g, t, r)r475] (open) 
and the mixing term (filled) in ea. (|66|) . Bottom: tr \^%{q,t,T)T (open) and the mixing term 
(filled) in eo.(|H7|). 
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. 20: Same qunatities for the proton, as in the previous figure 
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FIG. 21: The parity odd form factor F3 for neutron (top) and proton (bottom). Circles represent 
the result with 1^475 projection of ea. ()66|) . while triangles represent the one with T^^jj projection 
averaged over j = 1, 2, 3 in ea.(|67|). 
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FIG. 22: Effective mass plot for the parity-odd nucleon for \p\ = (top) and \p\ = 1 (bottom). 
Horizontal lines denote the central value (solid) and the errors (dotted) of the global fit. 
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